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Abstract  Electromagnetic  induction  (EMI:  10's  of  Hz  to  100's  of  kHz)  is  the  leading  technology  for 
discrimination  of  subsurface  metallic  targets  such  as  unexploded  ordnance  (UXO).  The  cleanup  problem 
requires  solution  of  remote  sensing  inverse  problem  inevitably  based  on  some  very  fast  forward 
algorithms  for  calculating  EMI.  The  forward  model  must  determine  responses  of  arbitrarily  complicated 
metallic  objects.  Here  a  very  fast  and  complete  forward  solution  system  is  presented,  based  on 
fundamental  mode  excitations.  For  a  given  target  (or  a  set  of  targets),  the  EMI  responses  to  fundamental 
modes  are  obtained  from  training  data  and  saved.  Any  realistic  excitation  field  is  then  decomposed  into  a 
limited  number  of  constituent  fundamental  modes  and  the  scatterer's  EMI  response  is  obtained  by 
superposition  of  the  fundamental  mode  solutions.  In  this  paper  we  define  the  fundamental  excitations 
explicitly  and  consider  their  rationale;  show  how  to  construct  any  particular  solutions  from  solutions  to 
the  fundamental  excitations;  and  focus  particularly  on  how  to  obtain,  retain,  express  responses  to  the 
fundamental  solutions  in  the  face  of  inherent  ill-conditioning. 

1  Introduction 


Cleanup  of  UXO  sites  is  extremely  expensive  because  of  high  false  alarm  rates.  Given  that  there  are  over 
10,000,000  acres  of  land  with  potential  UXO  hazard  in  the  US  alone,  reliable  techniques  for  UXO 
discrimination  are  urgently  needed.  As  an  inverse  problem,  UXO  discrimination  requires  a  fast  forward 
model,  i.e.  a  model  calculating  EMI  responses  for  prospective  targets  that  may  be  present.  Many 
numerical  techniques  are  available  for  EMI  calculation  in  the  magneto-quasistatic  EMI  regime,  but  most 
of  them  are  not  fast  enough  for  inversion  or  complete  enough  to  deal  with  realistic  objects.  One  of  the 
most  successful  models  is  the  dipole  model  [1,2],  in  which  one  approximates  a  target's  response  with  one 
or  a  number  of  infinitesimal  magnetic  dipoles,  each  responding  independently  to  the  local  value  of  the 
transmitted  ("primary")  field  that  stakes  it.  The  dipole  model  is  a  good  approximation  only  if  the 
observation  position  is  far  enough  from  the  target  and  if  the  parts  of  the  object  do  not  interact  significantly. 
However,  in  UXO  detection  and  discrimination  the  sensor  is  often  close  to  the  target,  and  we  have  shown 
elsewhere  that  interaction  effects  can  be  very  significant  [3,  4],  Fast  analytical  solution  for  a  spheroid  in 
EMI  frequencies  range  has  been  developed  recently  [5,6]  and  a  spheroid  model  was  shown  to  be  capable 
of  representing  effectively  some  geometrically  complex  objects  [7,8],  The  spheroid  model  has  two 
advantages:  (1)  It  takes  the  near  field  effect  into  account,  so  it  can  also  be  used  when  the  sensor  is  close; 
and  (2)  the  target  information  (geometrical  dimension,  orientation,  conductivity,  magnetic  permeability, 
etc.)  can  be  estimated  directly  from  model  parameters. 


For  more  complicated,  materially  heterogeneous  objects,  the  dipole  or  spheroid  models  are  not  wholly 
adequate,  and  detailed  numerical  solution  by  established  method  is  too  slow.  The  alternative  is  a 
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fundamental  mode  approach  [3,9,10].  We  choose  a  set  of  fundamental  excitation  modes  whose  linear 
combination  can  represent  an  arbitrary  excitation  field.  The  EMI  response  of  the  target  to  each 
fundamental  mode  is  obtained  and  stored.  Since  the  system  we  are  studying  is  linear,  if  some  general 
excitation  field  constitutes  a  particular  superposition  (linear  combination)  of  inputs,  the  response  will  be  a 
corresponding  superposition  of  outputs.  So  one  just  needs  to  build  a  library  for  known  UXOs,  where  the 
fundamental  solutions  (i.e.  the  response  of  each  UXO  to  fundamental  excitation  modes)  are  stored.  The 
EMI  response  of  any  UXO  to  an  arbitrary  excitation  can  then  be  easily  constructed  from  the  library  data. 
There  are  two  main  approaches  to  obtain  the  fundamental  solutions:  (1)  for  relatively  simple  object  whose 
components  (geometry  and  electromagnetic  properties  of  each  parts)  are  known,  the  fundamental  solution 
can  be  calculated  directly  through  numerical  simulation  [3].  (2)  For  the  more  general  case,  the 
fundamental  solution  can  be  obtained  indirectly  from  properly  designed  measurements  (so  called  training 
data).  A  major  difficulty  in  the  fundamental  mode  approach  is  that  the  excitation  fields  from  realistic  EMI 
sensors  are  usually  a  combination  of  the  fundamental  modes,  so  the  fundamental  solution  has  to  be  solved 
through  an  inversion  procedure,  which  often  suffers  from  ill-conditioning.  In  this  paper,  we  will  first 
introduce  a  selected  set  of  fundamental  modes  and  two  approaches  for  inferring  fundamental  solutions 
from  the  training  data.  Then  we  show  two  approaches  to  treat  the  ill-conditioning  problem  and 
demonstrate  some  results. 


2  Fundamental  mode  approach 

2.1  Expressing  arbitrary  primary  fields  using  fundamental  modes 

As  elaborated  elsewhere  [11,12],  in  the  EMI  band  we  can  usually  consider  the  magnetic  fields  in  the  air 
and  soil  to  be  irrotational,  with  the  result  that  they  can  be  expressed  as  gradients  of  a  scalar  potential  that 
satisfies  the  Laplace  equation.  Because  UXO  are  typically  elongated  bodies  of  revolution,  we  choose 
Associated  Legendre  functions  in  prolate  spheroidal  coordinates  as  the  fundamental  solution  forms.  On  a 
spheroid  surface  away  from  the  sensor,  the  primary  EMI  potential  field  can  be  written  as  a  summation  of 
fundamental  modes: 

v"  =  =  tw  (1) 

^  m  n—m  p— 0  k= 1 

where  here  and  in  what  follows  the  coordinate  system  is  centered  on  the  scatterer,  (tj,^,(|))  are  the 
standard  spheroidal  coordinates,  d  is  the  inter-focal  distance,  P"'  is  Associated  Legendre 

functions  of  the  first,  of  order  m  and  degree  n.  rpm((]))  is  cos (m</>)  for  p  =  0  and  is  sin (m</>)  for  p  = 
1.  The  coefficients  b  can  be  obtained  by  using  orthogonality  relations  for  the  Legendre  functions. 
For  the  purpose  of  simplicity,  we  rewrite  the  series  in  the  middle  of  (1)  using  a  condensed  subscript 
notation  on  the  right,  with  k  =  n-m  + 1  +  m  ■  nnax  +  p  ■  (mmax  + 1)  •  nmax  .  In  practice,  the  series  will  be 
truncated,  i.e.  we  will  only  consider  terms:  0<  p  <  1,  0  <  m  < mmax , m  <n<  m  +  nmax  -l.To  proceed  for  an 
arbitrary  object,  we  define  a  spheroidal  surface  S  surrounding  the  object,  more  or  less  conforming  to  its 
general  shape,  i.e.  we  choose  some  suitable  d  and  q0,  though  the  details  of  that  choice  are  not  critical. 

2.2  Constructing  general  responses  using  fundamental  solutions 


As  for  primary  field,  the  scattered  field  outside  of  the  target  (as  well  as  the  exterior  of  the  auxiliary 
spheroid)  is  irrotational  and  the  scattered  field  potential  function  satisfies  Laplace  equation.  For  each 

fundamental  excitation  ,  we  represent  the  corresponding  scattered  field  with  a  number  of  magnetic 
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charges  located  on  (or  inside)  the  auxiliary  spheroid  surface,  i.e.  the  scattered  potential  fields 
corresponding  to  excitation  y/j’1  at  location  R  are 


N  source 


MR)  =  Z 


,-=i  4/r//0|R-R,| 
and  the  total  scattered  magnetic  field  will  be 

P\ 


H(R)  =  YSL 


47Tju0  R-R, 


(RR) 


(2) 


(3) 


The  magnetic  charges  p[  are  obtained  by  the  methods  described  below  and  then  stored  before  attacking 
any  particular  application.  When  arbitrary  excitation  fields  are  decomposed  into  basic  modes  (i.e.  the  bk 
are  determined),  then  the  scattered  field  can  be  easily  computed  by  superposition  (i.e.  by  (3)  ).  One  of 
the  essential  steps  is  to  obtain  the  fundamental  solutions  in  the  first  place,  which  will  be  demonstrated  in 
the  next  section. 


2.3  Inferring  fundamental  solutions  from  training  data 

In  numerical  simulation,  EMI  response  of  a  specific  target  to  each  fundamental  excitation  y/^  can  be 

calculated  individually.  However,  it  is  not  the  case  in  real  measurement  because  our  instruments  do  not 
transmit  individual  fundamental  excitation  modes,  one  by  one.  The  field  from  realistic  sensor  is  usually  a 
combination  of  a  number  of  fundament  modes.  To  approach  this,  we  consider  below  the  particular  GEM- 
3  sensor  developed  by  Geophex  Ltd.  [13],  which  has  distinguished  itself  in  UXO  discrimination 
application. 


2.3.1  Simple  least  squares  approach 

The  GEM-3  we  are  currently  using  is  mono-static  and  measures  only  the  component  of  the  scattered  field 
normal  to  the  sensor  head,  so  the  forward  model  is  expressed  as: 


H 


4 Wo  R  -R, 


-(R.-R.J-Nj 


(4) 


Where  j  is  the  index  of  the  position  Rj  where  that  normal  field  was  measured,  and  N  j  is  the  direction 
normal  to  sensor  head  at  point  Rj . 


For  given  a  set  of  training  data  (well  designed  measurements  for  a  specific  target),  p'  can  be 
determined  by  minimizing 


msq=y.h,-h1  2=Z 


I 


4*9*0 


R-R 


(Rj-R,)  n,-h; 


(5) 


Using  a  least  square  error  approach,  we  can  build  the  normal  equations 
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A  q  =  b 

dMSQ  _ 

by  setting  — - =  U ,  i.e. 


Ai\i  ~  X 


(r,-r)n, 


}  y  k 


4^o  |R  ,  -  R, 

7  x-'  d  /'  (R-  R  )  N  ^ 

b,,  =e  «;  */■ 1 -  <w  = pi  s, 


bi  g  s 


R  -R 


R  -R 


Here  we  define 

£  _  I1  ^  1  =  n  ~ 111  + 1  +  m  ■  «max  +  f(»W  +  l)«max  +  2(l  -  l)(mmax  +  l)«max 

^  |o  ot/^gr 

The  magnetic  charge  magnitude  p[  (orq )  can  be  obtained  by  solving  (6). 


(6) 


(7) 


(8) 


2.3.2  Levenberg-Marquardt  (LM)  approach 

As  we  will  show  in  the  results  section,  the  inverse  problem  is  usually  highly  ill  conditioned.  For  this 
reason,  we  employ  LM  approach  [14],  where  the  magnetic  charge  magnitudes  qare  obtained  through 

iteration  with  a  regularization  factor  A, : 


qc*+D  =  qm  +  ( jTj  +  Xly  ^jT  JHd  _H) 


(9) 


Where  Hd  is  the  measurement  data,  H=H(q)  is  the  modeled  value,  and  the  Jacobian  matrix  J  = 


rJH 

3q 


3  Results 


3.1  Inherent  ill-conditioning 

In  general,  the  excitation  field  from  a  realistic  EMI  sensor  contains  mainly  a  small  number  of  basic  modes, 
with  the  contribution  of  other  modes  (usually  high  order  modes)  negligibly  small.  As  an  example,  Figure 

la  shows  the  coefficients  of  basic  modes  b  from  GEM-3  for  a  spheroid  with  equation 

(.r/0.03)-  +(z/0.18)~  =1  .  The  index  is  lined  up  according  to  the  formulation 

k  =  n- m  + 1  +  m ■  nmax  +  p ■  (rnmax  + 1)  ■  nmax  with  nmax  =  5,  mmax  =  3 .  The  first  term  bm) does  not  induce 
any  scattering.  So  in  this  case  the  dominant  terms  are  bom  (k=2,  corresponding  to  uniform  axial 
excitation),  bou  (k=6,  uniform  transverse  excitation)  and  several  other  small  terms.  All  the  other  terms  are 
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essentially  zero.  Figure  lb  is  similar  except  that  the  sensor  moves  along  a  cubic  box  surrounding  the 
spheroid  and  the  summation  of  the  coefficients  b  from  all  positions  are  plotted.  For  both  cases  we  see 
that  the  coefficients  of  many  terms  are  very  small  compared  to  the  dominant  ones. 


(a) 


(b) 


Figure  1  a)  coefficient  of  fundamental  modes  b  ,  when  the  spheroid  is  30  cm  from  the  40  cm  diameter  GEM-3 

sensor  and  is  titled  from  the  horizontal  at  45°  .  b)  Summation  of  | bpmn  |  over  all  observations,  as  the  sensor  moves 
over  a  90  cm  square  cube  surrounding  the  spheroid. 


3.2  Least  square  approach 


As  shown  in  previous  section,  the  training  data  contains  very  little  information  on  high  order  excitation 
modes.  The  consequence  is  that  EMI  data  in  this  set  of  training  data  is  not  sensitive  to  the  fundamental 
solution  of  these  terms,  i.e.  dH/dqk  is  close  to  zero,  so  that  the  matrix  A  in  (6)  is  singular  and  the 

equation  cannot  be  solved  accurately.  Fortunately,  these  high  order  terms  also  typically  contribute  very 
little  in  actual  applications,  so  that  the  accurate  fundamental  solution  of  these  terms  is  not  necessary.  This 
allows  us  to  simply  truncate  these  non-essential  terms  to  reduce  singularity.  In  this  least  squares  approach. 


we  define  the  criterion  of  truncation  as 


max  \b. 


i<  £  ,  where  max | bk |  is  the  maximum  of  all 


coefficients  except  k  —  1 .  Figure  2  shows  an  example  of  the  field  reconstruction  for  different  values  of  £. 


The  target  is  an  ellipsoid  given  by  (v/0.0l)" +(y/0.03)  +(z/0.18)~  =1  and  the  surrounding  spheroid 


is  (a/0.03)~  +(z/o.18)~  =  1  .  The  training  data  were  measured  along  a  cubic  box  with  surfaces 
at  x  =  ±0.45m  ,  y  =  ±0.45m  ,  and  z  =  +0.45m  ,  with  100  measurements  on  each  surface.  Instead  of 
actual  measurements,  the  training  data  in  this  paper  were  synthetic  ones  calculated  by  the  TSA  code  [1 1], 
To  better  understand  the  spheroid  position  relative  to  measurement  grid,  we  adopt  the  notation  of  Euler 
Angles  ( 0 ,  </>,ys)  [15].  This  box  measurement  is  equivalent  to  fixing  the  rectangular  grid  corresponding 

to  a  side  of  the  box  and  changing  the  center  location  (x0 ,  y0 ,  z0 )  and  orientation  of  the  target  and 
spheroid  as  follows:  x0  =  0  ,  y0  =  0  ,  z0  =  -45cm  ,  (01,</>1,l/fl)  =  ( 0,0,0),  (d2,02,l/f2)  =  (7U, 0,0), 
(6>3,03,^3)  =  (;r/2,0,0)  ,  (04,</>4,i//4)  =  (jt  I2,0,7T  IT)  ,  (05,</>5,y/5)  =  Or/2,0, ;r)  ,  (06,(/)6,y6)  = 

Or/  2,0,3#/  2) . 
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Figure  2  EMI  response  constructed  from  fundamental  solutions  (green  markers)  compared  to  synthetic  data  (circles) 


To  produce  the  fundamental  mode  solutions,  1 1  point  charges  are  distributed  uniformly  in  space  on  the 
suiTounding  spheroid  surface.  The  fundamental  solutions  were  obtained  from  the  training  data  via  (6).  To 
check  the  accuracy  of  the  inferred  fundamental  solutions,  we  calculated  the  EMI  response  along  rows  of 
observation  points  over  the  grid,  with  the  central  rows  passing  over  the  target,  producing  the  serrated 
response  pattern  in  Figure  2,  which  shows  results  for  the  phase  quadrature  component  of  47970  Hz 
excitation.  The  target  parameters  are  ( 6 ,  </),  y/)  =  (7t  1 2, 0, 0) ,  xf)  =  y0  =  0 ,  z0  =  —60 cm  .  The  truncation 

criteria  factors  are  £=  10~3,  10  s,  and  10~8  .  Results  show  that  we  can  obtain  reasonably  accurate 
fundamental  solutions  by  choosing  a  proper  truncation  criterion  (i.e.  10  3  and  10  5  in  this  case).  Our 
neglecting  the  truncated  terms  causes  the  remaining  small  error  in  the  best  results.  As  the  truncation 
criterion  becomes  too  small,  the  high  order  terms  begin  to  induce  singularity  and  may  destroy  the  whole 
solution. 


real  part,  f=330Hz 


Imaginary  part,  f=330Hz 


Figure  3  EMI  response  constructed  from  fundamental  solutions  (green  markers)  compared  to  measurements  (circles) 


Figure  3  showed  results  of  a  standard  81cm  UXO  for  both  real  and  imaginary  part  at  330Hz.  The 
UXO  is  centered  on  (-5,  -5,  -30)  cm  in  the  measurement  grid  coordinate,  with  a  horizontal 
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orientation:  (0,</),y/)  =  {n I2,n  12,0) .  The  training  data  are  3  set  of  measurements:  (1)  centerd  at  (- 
5,0.66,-39.48)cm,  angle  (#,^,yF)  =  (Jn  120, n  12,0)  (2)  (-5,0.66,-39.48)cm,  angle  ( 0,(f),y/)  = 

{29  k  I  26,2k  1 2,0)  (3)  center  at  (-5,  -5,-20)  cm,  angle  {0,(f),y/)  =  {k  ! 2,k/2, 0) .  The  fundamental 

solutions  are  calculated  with  truncation  parameter  £  =  1CT4  and  the  modeled  data  match  with  measured 
ones  pretty  well.  It  should  be  noticed  that  the  testing  data  is  close  to  one  set  of  training  data  (set  3).  For 
more  complete  test  of  the  model  accuracy,  more  general  measurements  need  to  be  done  for  comparison. 

3.3  LM  approach 

Instead  of  simply  truncating  the  “non-essential”  terms  (which  is  equivalent  to  setting  their  solution  to  be 
zero),  in  the  LM  approach  we  keep  these  terms  and  the  iterative  procedure  will  assign  them  some 
reasonable  values,  which  produces  more  accurate  solutions  than  simple  truncation.  This  regularization 
factor  X  makes  sure  that  their  value  won’t  become  so  erratic  as  to  cause  instability.  For  the  same  set  of 
training  data  as  in  section  3.2,  we  set  an  initial  value  for  X,  reducing  it  as  the  iteration  proceeds.  Once  the 
error  term  stops  decreasing,  we  set  the  regularization  factor  back  to  initial  value  and  iterate  again  [16]. 
Experience  shows  that  the  solution  is  not  very  sensitive  to  A,r  and  results  are  more  accurate  and  stable  than 
those  from  the  least  squares  approach.  Some  examples  are  shown  in  as  shown  in  Figure  4  . 


0  20  40  60  80  100  0  20  40  60  80  100 

Ns.  of  observations  Na  of  observations 


Na  of  observations 


Figure  4.  EMI  response  constmcted  from  fundamental  solutions  (green  markers)  compared  to  synthetic  data 
(circles),  based  on  the  LM  approach. 


4  Concluding  discussion 

A  fundamental  mode  approach  was  developed  for  fast  calculation  of  the  complete  EMI  response  of 
complicated  3D  objects.  For  prospective  targets,  the  fundamental  solutions  are  obtained  and  stored  in 
library.  The  response  of  these  targets  to  arbitrary  EMI  excitation  can  then  be  constructed  quickly  by 
simple  superposition  of  the  fundamental  solutions.  Two  approaches  to  infer  fundamental  solutions  from 
training  data  (usually  measurement  data)  were  introduced:  simple  least  squares  error  minimization  and  a 
Levenberg-Marquardt  approach.  Both  approach  give  reasonable  results  if  inversion  control  parameters  are 
properly  used.  While  the  preliminary  results  are  promising,  further  study  is  needed  to  make  the  approach 
more  efficient  and  practical.  More  realistic  targets  need  to  be  studied  to  determine  the  criteria  for 
truncation  and  regularization  parameters.  We  also  seek  better  ways  of  distributing  the  magnetic  charges  to 
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reduce  the  number  of  unknowns  and  the  computational  cost.  A  better  design  of  training  data  (especially 
bi-static  data  or  data  for  three  components  of  H)  may  also  help  to  reduce  ill-conditioning  problems. 
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